/*************************************************************************
 *
 * This file is part of the SAMRAI distribution.  For full copyright
 * information, see COPYRIGHT and LICENSE.
 *
 * Copyright:     (c) 1997-2020 Lawrence Livermore National Security, LLC
 * Description:   Utilitiy for breaking boxes during partitioning.
 *
 ************************************************************************/

#ifndef included_mesh_BalanceBoxBreaker_C
#define included_mesh_BalanceBoxBreaker_C

#include "SAMRAI/mesh/BalanceBoxBreaker.h"

#include "SAMRAI/hier/BoxUtilities.h"
#include "SAMRAI/mesh/BalanceUtilities.h"
#include "SAMRAI/tbox/Utilities.h"
#include "SAMRAI/tbox/TimerManager.h"

#if !defined(__BGL_FAMILY__) && defined(__xlC__)
/*
 * Suppress XLC warnings
 */
#pragma report(disable, CPPC5334)
#pragma report(disable, CPPC5328)
#endif

namespace SAMRAI {
namespace mesh {

// Round a to the nearest higher integer divisible by b.  This should work even for a < 0.
#define ROUND_TO_HI(a, b) ((a) - ((((a) % (b)) - (b)) % (b)))
// Round a to the nearest lower integer divisible by b.  This should work even for a < 0.
#define ROUND_TO_LO(a, b) ((a) - ((((a) % (b)) + (b)) % (b)))

BalanceBoxBreaker::BalanceBoxBreaker(
   const PartitioningParams& pparams,
   bool print_break_steps):
   d_pparams(&pparams),
   d_print_break_steps(print_break_steps)
{
   setTimers();
}

BalanceBoxBreaker::BalanceBoxBreaker(
   const BalanceBoxBreaker& other):
   d_pparams(other.d_pparams),
   d_print_break_steps(other.d_print_break_steps),
   t_break_off_load(other.t_break_off_load),
   t_find_bad_cuts(other.t_find_bad_cuts)
{
}

/*
 *************************************************************************
 * Master method for breaking off a load.
 *
 * Try different heuristics and pick the "best" way to break off a
 * load.  The best is defined as the one with the lowest combined
 * penalty.
 *
 * This method always return a breakage if at all possible, without
 * considering whether the break should be used.  For example,
 * requesting breakage of 1 cell in a 100x100 box might return a
 * breakage of a 100-cells sliver!
 *
 * Return whether a successful break was made.
 *************************************************************************
 */
bool
BalanceBoxBreaker::breakOffLoad(
   hier::BoxContainer& breakoff,
   hier::BoxContainer& leftover,
   double& brk_load,
   const hier::Box& box,
   double box_load,
   const std::vector<double>& corner_weights,
   double ideal_load_to_break,
   double low_load,
   double high_load,
   double threshold_width) const
{
   TBOX_ASSERT(ideal_load_to_break > 0);

   /*
    * NOTE: We need in this method a way to weigh the
    * value of proximity to the ideal breakoff vs the
    * increased area of the cuts.  However, the weight
    * given to area-optimized cuts should be considered
    * only with real application performance data.
    *
    * NOTE: We can compute the amount of new box boundaries
    * generated by computing the box boundary before and
    * after, and subtracting.  Easier than reconstructing
    * the cuts from the box definitions.
    *
    * NOTE: We should weight negatively the production of
    * high surface-to-volume boxes.
    */

   t_break_off_load->start();

   breakoff.clear();
   leftover.clear();

   /*
    * To avoid repeated computations of bad cuts,
    * we precompute bad_cuts here to provide to
    * methods that actually use the information.
    */
   std::vector<std::vector<bool> > bad_cuts(box.getDim().getValue());
   t_find_bad_cuts->start();
   hier::BoxUtilities::findBadCutPoints(
      bad_cuts,
      box,
      d_pparams->getDomainBoxes(box.getBlockId()),
      d_pparams->getBadInterval());
   t_find_bad_cuts->stop();

   // Penalty for not transfering ideal load.
   brk_load = 0;
   double best_balance_penalty = computeBalancePenalty(ideal_load_to_break);
   double best_width_score = computeWidthScore(box.numberCells(), threshold_width);

   if (d_print_break_steps) {
      tbox::plog.unsetf(std::ios::fixed | std::ios::scientific);
      tbox::plog.precision(6);
      tbox::plog << "      breakOffLoad: pre-break imbalance: " << ideal_load_to_break
                 << " balance penalty: " << best_balance_penalty
                 << " width score: " << best_width_score
                 << " threshold_width: " << threshold_width
                 << std::endl;
   }

   TrialBreak best_trial(*d_pparams, threshold_width,
                         box, box_load, corner_weights, bad_cuts,
                         ideal_load_to_break, low_load, high_load);

   /*
    * Try planar break.
    */

   TrialBreak planar_trial(*d_pparams, threshold_width,
                           box, box_load, corner_weights, bad_cuts,
                           ideal_load_to_break, low_load, high_load);

   if (breakOffLoad_planar(planar_trial)) {

      if (d_print_break_steps) {
         // Should move this if-block into breakOffLoad_planar.
         tbox::plog.unsetf(std::ios::fixed | std::ios::scientific);
         tbox::plog.precision(6);
         tbox::plog << "      breakOffLoad_planar broke off "
                    << planar_trial.d_breakoff_load << " / " << ideal_load_to_break
                    << " from " << box << '|'
                    << box.numberCells() << '|'
                    << box.size() << " into "
                    << planar_trial.d_breakoff.size()
                    << " breakoff: ";
         for (hier::BoxContainer::const_iterator bi =
                 planar_trial.d_breakoff.begin();
              bi != planar_trial.d_breakoff.end();
              ++bi) {
            tbox::plog << " " << *bi << '|' << bi->numberCells() << '|'
                       << bi->size();
         }
         tbox::plog << "\n        and " << planar_trial.d_leftover.size()
                    << " leftover boxes:";
         for (hier::BoxContainer::const_iterator bi =
                 planar_trial.d_leftover.begin();
              bi != planar_trial.d_leftover.end();
              ++bi) {
            tbox::plog << " " << *bi << '|' << bi->numberCells() << '|'
                       << bi->size();
         }
         tbox::plog << "\n        imbalance: "
                    << (planar_trial.d_breakoff_load - ideal_load_to_break)
                    << planar_trial
                    << std::endl;
      }

      if (planar_trial.improvesOver(best_trial)) {
         best_trial.swap(planar_trial);
         if (d_print_break_steps) {
            tbox::plog << " breakOffLoad: chose the above trial break" << std::endl;
         }
      }

   } else {
      if (d_print_break_steps) {
         tbox::plog << "      breakOffLoad_planar could not break "
                    << ideal_load_to_break << " from " << box
                    << '/' << box.numberCells()
                    << '/' << box.numberCells().getProduct()
                    << std::endl;
      }
   }

   /*
    * Try cubic break.
    */

   TrialBreak cubic_trial(*d_pparams, threshold_width,
                          box, box_load, corner_weights, bad_cuts,
                          ideal_load_to_break, low_load, high_load);

   if (breakOffLoad_cubic(cubic_trial)) {

      if (d_print_break_steps) {
         // Should move this if-block into breakOffLoad_cubic.
         tbox::plog.unsetf(std::ios::fixed | std::ios::scientific);
         tbox::plog.precision(6);
         tbox::plog << "      breakOffLoad_cubic broke off "
                    << cubic_trial.d_breakoff_load << " / " << ideal_load_to_break
                    << " from " << box << '|'
                    << box.numberCells() << '|'
                    << box.size() << " into "
                    << cubic_trial.d_breakoff.size()
                    << " breakoff: ";
         for (hier::BoxContainer::const_iterator bi =
                 cubic_trial.d_breakoff.begin();
              bi != cubic_trial.d_breakoff.end();
              ++bi) {
            tbox::plog << " " << *bi << '|' << bi->numberCells() << '|'
                       << bi->size();
         }
         tbox::plog << "\n        and " << cubic_trial.d_leftover.size()
                    << " leftover boxes:";
         for (hier::BoxContainer::const_iterator bi =
                 cubic_trial.d_leftover.begin();
              bi != cubic_trial.d_leftover.end();
              ++bi) {
            tbox::plog << " " << *bi << '|' << bi->numberCells() << '|'
                       << bi->size();
         }
         tbox::plog << "\n        imbalance: "
                    << (cubic_trial.d_breakoff_load - ideal_load_to_break)
                    << cubic_trial
                    << std::endl;
      }

      if (cubic_trial.improvesOver(best_trial)) {
         best_trial.swap(cubic_trial);
         if (d_print_break_steps) {
            tbox::plog << " breakOffLoad: chose the above trial break" << std::endl;
         }
      } else {
         if (d_print_break_steps) {
            tbox::plog << " breakOffLoad: cubic break did not improve:"
                       << "\n   cubic_trial:" << cubic_trial
                       << "\n   best_trial:" << best_trial
                       << std::endl;
         }
      }

   } else {
      if (d_print_break_steps) {
         tbox::plog << "      breakOffLoad_cubic could not break "
                    << ideal_load_to_break << " from " << box
                    << '/' << box.numberCells()
                    << '/' << box.numberCells().getProduct()
                    << std::endl;
      }
   }

   if (d_print_break_steps) {
      if (!best_trial.d_breakoff.empty()) {
         tbox::plog << "      breakOffLoad: broke:"
                    << best_trial << std::flush;
      } else {
         tbox::plog << "      breakOffLoad: no break!"
                    << std::endl;
      }
   }

   if (!best_trial.d_breakoff.empty()) {
      breakoff.swap(best_trial.d_breakoff);
      leftover.swap(best_trial.d_leftover);
      brk_load = best_trial.d_breakoff_load;
   }

   t_break_off_load->stop();
   return !breakoff.empty();
}

/*
 *************************************************************************
 * Break off a load from a box by making a single planar cut across
 * the box's longest direction.
 *
 * Currently assuming uniform load of one unit per cell.
 *
 * Return whether a successful break was made.
 *************************************************************************
 */
bool
BalanceBoxBreaker::breakOffLoad_planar(TrialBreak& trial) const
{
   const tbox::Dimension dim(trial.d_whole_box.getDim());

   if (d_print_break_steps) {
      tbox::plog << "        breakOffLoad_planar attempting to break "
                 << trial.d_ideal_load << " from Box "
                 << trial.d_whole_box << trial.d_whole_box.numberCells() << '|'
                 << trial.d_whole_box.size()
                 << " min_size=" << d_pparams->getMinBoxSize() << std::endl;
   }

   trial.d_breakoff.clear();
   trial.d_leftover.clear();

   if (trial.d_whole_box.size() <= d_pparams->getMinimumCellRequest()) {
      return false;
   }

   const hier::IntVector box_dims = trial.d_whole_box.numberCells();

   if (trial.d_whole_box_load <= trial.d_ideal_load) {
      // Easy: break off everything.
      trial.computeBreakData(trial.d_whole_box);
      if (d_print_break_steps) {
         tbox::plog << "      breakOffLoad_planar broke off entire Box "
                    << trial.d_whole_box
                    << std::endl;
      }
      return true;
   }

   /*
    * Determine ordering of box_dims from shortest to longest.
    */
   hier::IntVector sorted_dirs(dim);
   sorted_dirs.sortIntVector(box_dims);

   bool sufficient_brk_load = false;
   const hier::BlockId& block_id = trial.d_whole_box.getBlockId();
   hier::IntVector cut_factor(
      d_pparams->getCutFactor().getBlockVector(block_id));

   TrialBreak trial1(trial);

   for (tbox::Dimension::dir_t d1 = 1; d1 <= dim.getValue() && !sufficient_brk_load; ++d1) {
      const tbox::Dimension::dir_t d = static_cast<tbox::Dimension::dir_t>(dim.getValue() - d1);

      /*
       * Search directions from longest to shortest because we prefer
       * to break across longer directions.
       */
      const tbox::Dimension::dir_t brk_dir = static_cast<tbox::Dimension::dir_t>(sorted_dirs(d));

      const std::vector<bool>& bad = trial.d_bad_cuts[brk_dir];

      const double ideal_fraction = static_cast<double>(trial.d_ideal_load) /
                                    static_cast<double>(trial.d_whole_box_load);
      const double ideal_cut_length = ideal_fraction *
                                      static_cast<double>(box_dims(brk_dir));

      /*
       * Try 4 different cuts for direction brk_dir:
       * upper/lower: should we break off the upper end or lower end of the box.
       * hi/lo: should we round the break plane to the high or low side.
       *
       * plane refers to the index of the mesh line where the cut is.
       *
       * offset refers to a cut plane index relative to the box's lower corner.
       */

      // Ideal cut planes, not necessarily coincident with a grid line.
      const double ideal_upper_cut_offset = trial.d_whole_box.numberCells(brk_dir)
         - ideal_cut_length;
      const double ideal_lower_cut_offset = ideal_cut_length;

      // Compute valid cut planes on high and low sides of upper cut plane.
      int lo_upper_cut_plane = trial.d_whole_box.lower() (brk_dir) + int(ideal_upper_cut_offset);
      int hi_upper_cut_plane = trial.d_whole_box.lower() (brk_dir) + int(ideal_upper_cut_offset)
         + 1;
      lo_upper_cut_plane = ROUND_TO_LO(lo_upper_cut_plane, cut_factor(brk_dir));
      hi_upper_cut_plane = ROUND_TO_HI(hi_upper_cut_plane, cut_factor(brk_dir));
      while (lo_upper_cut_plane > trial.d_whole_box.lower() (brk_dir) &&
             bad[lo_upper_cut_plane - trial.d_whole_box.lower() (brk_dir)]) {
         lo_upper_cut_plane -= cut_factor(brk_dir);
      }
      while (hi_upper_cut_plane < trial.d_whole_box.upper() (brk_dir) + 1 &&
             bad[hi_upper_cut_plane - trial.d_whole_box.lower() (brk_dir)]) {
         hi_upper_cut_plane += cut_factor(brk_dir);
      }

      // Compute valid cut planes on high and low sides of lower cut plane.
      int lo_lower_cut_plane = trial.d_whole_box.lower() (brk_dir) + int(ideal_lower_cut_offset);
      int hi_lower_cut_plane = trial.d_whole_box.lower() (brk_dir) + int(ideal_lower_cut_offset)
         + 1;
      lo_lower_cut_plane = ROUND_TO_LO(lo_lower_cut_plane, cut_factor(brk_dir));
      hi_lower_cut_plane = ROUND_TO_HI(hi_lower_cut_plane, cut_factor(brk_dir));
      while (lo_lower_cut_plane > trial.d_whole_box.lower() (brk_dir) &&
             bad[lo_lower_cut_plane - trial.d_whole_box.lower() (brk_dir)]) {
         lo_lower_cut_plane -= cut_factor(brk_dir);
      }
      while (hi_lower_cut_plane < trial.d_whole_box.upper() (brk_dir) + 1 &&
             bad[hi_lower_cut_plane - trial.d_whole_box.lower() (brk_dir)]) {
         hi_lower_cut_plane += cut_factor(brk_dir);
      }

      if (d_print_break_steps) {
         tbox::plog << "breakOffLoad_planar considering cuts along " << d << ": "
                    << "  ideal_cut_length=" << ideal_cut_length
                    << "  ideal_upper_cut_offset=" << ideal_upper_cut_offset
                    << "  lo_lower_cut_plane=" << lo_lower_cut_plane
                    << "  hi_lower_cut_plane=" << hi_lower_cut_plane
                    << "  lo_upper_cut_plane=" << lo_upper_cut_plane
                    << "  hi_upper_cut_plane=" << hi_upper_cut_plane
                    << std::endl;
      }

      if (lo_lower_cut_plane - trial.d_whole_box.lower() (brk_dir) >
          d_pparams->getMinBoxSize() (brk_dir) &&
          trial.d_whole_box.upper() (brk_dir) + 1 - lo_lower_cut_plane >
          d_pparams->getMinBoxSize() (brk_dir)) {

         hier::Box brk_box(trial.d_whole_box);
         brk_box.setUpper(brk_dir, lo_lower_cut_plane - 1);
         trial1.computeBreakData(brk_box);
         if (trial1.improvesOver(trial)) {
            trial.swap(trial1);
            if (d_print_break_steps) {
               tbox::plog << "breakOffLoad_planar choosing dir " << brk_dir
                          << " lo_lower_cut_plane " << ":" << trial << std::endl;
            }
         }
      }

      if ((hi_lower_cut_plane - trial.d_whole_box.lower() (brk_dir) >
           d_pparams->getMinBoxSize() (brk_dir) &&
           trial.d_whole_box.upper() (brk_dir) + 1 - hi_lower_cut_plane >
           d_pparams->getMinBoxSize() (brk_dir)) ||
          hi_lower_cut_plane >= trial.d_whole_box.upper() (brk_dir) + 1) {

         hier::Box brk_box(trial.d_whole_box);
         brk_box.setUpper(brk_dir,
            tbox::MathUtilities<int>::Min(hi_lower_cut_plane - 1,
               trial1.d_whole_box.upper(brk_dir)));
         trial1.computeBreakData(brk_box);
         if (trial1.improvesOver(trial)) {
            trial.swap(trial1);
            if (d_print_break_steps) {
               tbox::plog << "breakOffLoad_planar choosing dir " << brk_dir
                          << " hi_lower_cut_plane " << ":" << trial << std::endl;
            }
         }
      }

      if ((trial.d_whole_box.upper() (brk_dir) + 1 - lo_upper_cut_plane >
           d_pparams->getMinBoxSize() (brk_dir) &&
           lo_upper_cut_plane - trial.d_whole_box.lower() (brk_dir) >
           d_pparams->getMinBoxSize() (brk_dir)) ||
          lo_upper_cut_plane <= trial.d_whole_box.lower() (brk_dir)) {

         hier::Box brk_box(trial.d_whole_box);
         brk_box.setLower(brk_dir,
            tbox::MathUtilities<int>::Max(lo_upper_cut_plane,
               trial1.d_whole_box.lower(brk_dir)));
         trial1.computeBreakData(brk_box);
         if (trial1.improvesOver(trial)) {
            trial.swap(trial1);
            if (d_print_break_steps) {
               tbox::plog << "breakOffLoad_planar choosing dir " << brk_dir
                          << " lo_upper_cut_plane " << ":" << trial << std::endl;
            }
         }
      }

      if (trial.d_whole_box.upper() (brk_dir) + 1 - hi_upper_cut_plane >
          d_pparams->getMinBoxSize() (brk_dir) &&
          hi_upper_cut_plane - trial.d_whole_box.lower() (brk_dir) >
          d_pparams->getMinBoxSize() (brk_dir)) {

         hier::Box brk_box(trial.d_whole_box);
         brk_box.setLower(brk_dir, hi_upper_cut_plane);
         trial1.computeBreakData(brk_box);
         if (trial1.improvesOver(trial)) {
            trial.swap(trial1);
            if (d_print_break_steps) {
               tbox::plog << "breakOffLoad_planar choosing dir " << brk_dir
                          << " hi_upper_cut_plane " << ":" << trial << std::endl;
            }
         }
      }

      sufficient_brk_load = (trial.d_breakoff_load >= trial.d_low_load) &&
         (trial.d_breakoff_load <= trial.d_high_load);

   } // d-loop

#ifdef DEBUG_CHECK_ASSERTIONS
   for (hier::BoxContainer::iterator bi = trial.d_breakoff.begin();
        bi != trial.d_breakoff.end();
        ++bi) {
      const hier::Box& b = *bi;
      const hier::IntVector s = b.numberCells();
      for (int d = 0; d < dim.getValue(); ++d) {
         if (((s(d) < d_pparams->getMinBoxSize() (d)) && (s(d) != box_dims(d))) ||
             (s(d) > box_dims(d))) {
            TBOX_ERROR("BalanceBoxBreaker library error:\n"
               << "breakoff box " << b << ", size " << s
               << "\nis not between the min size " << d_pparams->getMinBoxSize()
               << "\nand the original box size " << box_dims << "\n"
               << "break box size " << trial.d_breakoff.front().numberCells() << "\n"
               << "ideal brk load " << trial.d_ideal_load);
         }
      }
   }
   for (hier::BoxContainer::iterator bi = trial.d_leftover.begin();
        bi != trial.d_leftover.end();
        ++bi) {
      const hier::Box& b = *bi;
      const hier::IntVector s = b.numberCells();
      for (int d = 0; d < dim.getValue(); ++d) {
         if (((s(d) < d_pparams->getMinBoxSize() (d)) && (s(d) != box_dims(d))) ||
             (s(d) > box_dims(d))) {
            TBOX_ERROR("BalanceBoxBreaker library error:\n"
               << "leftover box " << b << ", size " << s
               << "\nis not between the min size " << d_pparams->getMinBoxSize()
               << "\nand the original box size " << box_dims << "\n"
               << "break box size " << trial.d_breakoff.front().numberCells() << "\n"
               << "ideal brk load " << trial.d_ideal_load);
         }
      }
   }
#endif
   if (d_print_break_steps) {
      tbox::plog << "        breakOffLoad_planar returning." << std::endl;
   }

   return !trial.d_breakoff.empty();
}

/*
 *************************************************************************
 * Break off a load from a box by making a box cut that is as close
 * to cubic as feasible.
 *
 * Currently assuming uniform load of one unit per cell.
 *
 * Return whether a successful break was made.
 *
 * Differs from breakOffLoad in that it will always
 * performs a break and if needed, break off more than
 * the ideal.  The calling method should take this into
 * account.
 *************************************************************************
 */
bool
BalanceBoxBreaker::breakOffLoad_cubic(TrialBreak& trial) const
{
   if (trial.d_whole_box.size() <= d_pparams->getMinimumCellRequest()) {
      return false;
   }

   const tbox::Dimension dim(trial.d_whole_box.getDim());

   const hier::IntVector box_dims(trial.d_whole_box.numberCells());

   if (trial.d_ideal_load >= trial.d_whole_box_load) {
      // Easy: break off everything.
      trial.computeBreakData(trial.d_whole_box);
      if (d_print_break_steps) {
         tbox::plog << "      breakOffLoad_cubic broke off entire Box "
                    << trial.d_whole_box
                    << std::endl;
      }
      return true;
   }

   if (trial.d_ideal_load > 0.5 * trial.d_whole_box_load &&
       !tbox::MathUtilities<double>::equalEps(trial.d_ideal_load, 0.5 * trial.d_whole_box_load)) {
      /*
       * This algorithm is better when breaking off a small portion.
       * Since the ideal is a bigger portion, switch breakoff with leftover.
       */
      if (d_print_break_steps) {
         tbox::plog
         << "      breakOffLoad_cubic reversing direction to break "
         << (static_cast<double>(box_dims.getProduct()) - trial.d_ideal_load)
         << " instead of " << trial.d_ideal_load << " / "
         << box_dims.getProduct() << std::endl;
      }
      TrialBreak reversed_trial(trial, true);
      bool success =
         breakOffLoad_cubic(reversed_trial);
      if (success) {
         trial.swapWithReversedTrial(reversed_trial);
      }
      return success;
   }

   if (d_print_break_steps) {
      tbox::plog << "      breakOffLoad_cubic attempting to break "
                 << trial.d_ideal_load << " from Box "
                 << trial.d_whole_box << trial.d_whole_box.numberCells() << '|'
                 << trial.d_whole_box.size()
                 << " min_size=" << d_pparams->getMinBoxSize() << std::endl;
   }

   trial.d_breakoff.clear();
   trial.d_leftover.clear();
   trial.d_breakoff_load = 0.0;

   const hier::IntVector& one_vec = hier::IntVector::getOne(dim);
   const hier::IntVector& zero_vec = hier::IntVector::getZero(dim);

   hier::Box best_breakoff_box(dim);
   hier::IntVector best_breakoff_size = zero_vec;
   double best_breakoff_load = 0;

   const hier::BlockId& block_id = trial.d_whole_box.getBlockId();
   const hier::IntVector& cut_factor(
      d_pparams->getCutFactor().getBlockVector(block_id));

   /*
    * We consider 2^dim boxes grown from the incoming box's corners.
    * Imagine 2 cutting planes per dimension, perpendicular to each
    * axis, an upper cut and a lower cut.  The box is cut into 3 parts
    * in each direction, so the 2*dim planes divide the box into 3^dim
    * boxes.
    *
    * 2D example:
    *
    *       +-----------------+
    *       |    |        |   |
    *       |----+--------+---|
    *       |    |        |   |
    *       |    |        |   |
    *       |    |        |   |
    *       |    |        |   |
    *       |----+--------+---|
    *       |    |        |   |
    *       +-----------------+
    *
    * upper_intersection is the point where all upper cuts intersect.
    * lower_intersection is the point where all lower cuts intersect.
    * We only consider the cut boxes that touch the incoming box's
    * corners.  Using the other boxes result in too much fragmentation
    * of the incoming box.
    */
   hier::IntVector brk_size(d_pparams->getMinBoxSize());
   brk_size.max(cut_factor);
   brk_size.min(box_dims);

   /*
    * Make sure brk_size is a multiple of cut_factor.
    */
   for (int d = 0; d < dim.getValue(); ++d) {
      if (brk_size(d) % cut_factor(d) != 0) {
         brk_size(d) = ((brk_size(d) / cut_factor(d)) + 1) * cut_factor(d);
      }
   }

   /*
    * Initialize the cut plane intersections.  We will grow the
    * corner boxes by gradually moving the intersections away
    * from their initial location.
    */
   hier::IntVector lower_intersection(trial.d_whole_box.lower() + d_pparams->getMinBoxSize());
   hier::IntVector upper_intersection(
      trial.d_whole_box.upper() - d_pparams->getMinBoxSize() + one_vec);
   for (int d = 0; d < dim.getValue(); ++d) {
      lower_intersection(d) = ROUND_TO_HI(lower_intersection(d),
            cut_factor(d));
      upper_intersection(d) = ROUND_TO_LO(upper_intersection(d),
            cut_factor(d));
   }

   const int num_corners = 1 << dim.getValue();

   for (int bn = 0; bn < num_corners; ++bn) {

      if (d_print_break_steps) {
         tbox::plog << "Examining corner box " << bn << std::endl;
      }

      /*
       * Compute initial box at corner bn and its expansion rate.
       */
      hier::Box corner_box(trial.d_whole_box);
      hier::IntVector corner_box_size = zero_vec;
      double corner_box_load = 0;
      hier::IntVector expansion_rate(dim);

      for (int d = 0; d < dim.getValue(); ++d) {

         // In direction d, does corner_box touch the upper (vs lower) side of box:
         int touches_upper_side = bn & (1 << d);

         if (touches_upper_side) {
            corner_box.setLower(static_cast<hier::Box::dir_t>(d),
               upper_intersection(d));
            if (corner_box.lower(static_cast<hier::Box::dir_t>(d))
                - trial.d_whole_box.lower(static_cast<hier::Box::dir_t>(d)) <
                d_pparams->getMinBoxSize() (d)) {
               corner_box.setLower(static_cast<hier::Box::dir_t>(d),
                  trial.d_whole_box.lower(static_cast<hier::Box::dir_t>(d)));
            }
            expansion_rate(d) = -cut_factor(d);
         }
         else {
            corner_box.setUpper(static_cast<hier::Box::dir_t>(d),
               lower_intersection(d) - 1);
            if ( trial.d_whole_box.upper()(d) - corner_box.upper()(d) < d_pparams->getMinBoxSize()(d) ) {
               corner_box.setUpper(static_cast<hier::Box::dir_t>(d),
                  trial.d_whole_box.upper(static_cast<hier::Box::dir_t>(d)));
            }
            expansion_rate(d) = cut_factor(d);
         }

      }

      corner_box_size = corner_box.numberCells();
      corner_box_load = trial.computeBreakOffLoad(corner_box);

      if (d_print_break_steps) {
         tbox::plog << "Initial corner box " << bn << " is " << corner_box
                    << corner_box.numberCells() << '|' << corner_box.size() << std::endl;
      }

      int break_acceptance_flags[4] = { 0, 0, 0, 0 };

      if (corner_box.size() >= d_pparams->getMinimumCellRequest() &&
          BalanceUtilities::compareLoads(
             break_acceptance_flags, best_breakoff_load,
             corner_box_load, trial.d_ideal_load,
             trial.d_low_load, trial.d_high_load, *d_pparams)) {
         best_breakoff_box = corner_box;
         best_breakoff_size = corner_box_size;
         best_breakoff_load = corner_box_load;
         if (d_print_break_steps) {
            tbox::plog << "best_breakoff_box is now box " << bn << " " << best_breakoff_box
                       << ", best_breakoff_size = " << best_breakoff_size
                       << ", best_breakoff_load = " << best_breakoff_load
                       << std::endl;
         }
         if (best_breakoff_load >= trial.d_ideal_load) {
            break;
         }
      }

      /*
       * growable: whether corner_box_size can be grown without
       * breaking off too much.
       */
      hier::IntVector growable(dim, 1);
      for (int d = 0; d < dim.getValue(); ++d) {
         growable[d] = corner_box_size[d] < box_dims[d];
      }

      while (corner_box_load < trial.d_ideal_load) {
         /*
          * The while loop gradually increases corner_box to bring
          * its size closer to trial.d_ideal_load.  Stop loop when
          * its size is in the acceptable range or if increasing
          * it only takes it farther from trial.d_ideal_load.
          *
          * Select inc_dir, the direction to expand corner_box.  Use the
          * smallest direction that is still allowed to grow.
          */

         int inc_dir = -1;
         for (int d = 0; d < dim.getValue(); ++d) {
            if (growable(d) &&
                (inc_dir == -1 || corner_box_size(d) < corner_box_size(inc_dir)))
               inc_dir = d;
         }
         if (inc_dir == -1) break;  // No growable direction.

         TBOX_ASSERT(corner_box_size(inc_dir) < box_dims(inc_dir));

         /*
          * Grow corner_box, but keep within boundary of box and
          * prevent remainder from violating min size.  Update growability.
          */
         if (expansion_rate(inc_dir) > 0) {
            corner_box.setUpper(static_cast<hier::Box::dir_t>(inc_dir),
               tbox::MathUtilities<int>::Min(
                  corner_box.upper(static_cast<hier::Box::dir_t>(inc_dir)) + expansion_rate(inc_dir),
                  trial.d_whole_box.upper(static_cast<hier::Box::dir_t>(inc_dir))));
            if (trial.d_whole_box.upper() (inc_dir) - corner_box.upper() (inc_dir) <
                d_pparams->getMinBoxSize() (inc_dir)) {
               corner_box.setUpper(static_cast<hier::Box::dir_t>(inc_dir),
                  trial.d_whole_box.upper(static_cast<hier::Box::dir_t>(inc_dir)));
            }
            growable(inc_dir) = corner_box.upper() (inc_dir) < trial.d_whole_box.upper() (inc_dir);
         } else {
            corner_box.setLower(static_cast<hier::Box::dir_t>(inc_dir),
               tbox::MathUtilities<int>::Max(
                  corner_box.lower(static_cast<hier::Box::dir_t>(inc_dir)) + expansion_rate(inc_dir),
                  trial.d_whole_box.lower(static_cast<hier::Box::dir_t>(inc_dir))));
            if (corner_box.lower() (inc_dir) - trial.d_whole_box.lower() (inc_dir) <
                d_pparams->getMinBoxSize() (inc_dir)) {
               corner_box.setLower(static_cast<hier::Box::dir_t>(inc_dir),
                  trial.d_whole_box.lower(static_cast<hier::Box::dir_t>(inc_dir)));
            }
            growable(inc_dir) = corner_box.lower() (inc_dir) > trial.d_whole_box.lower() (inc_dir);
         }

         corner_box_size = corner_box.numberCells();

         if (corner_box.size() < d_pparams->getMinimumCellRequest()) {
            continue;
         }

         corner_box_load = trial.computeBreakOffLoad(corner_box);

         const bool accept_break = BalanceUtilities::compareLoads(
               break_acceptance_flags, best_breakoff_load, corner_box_load,
               trial.d_ideal_load, trial.d_low_load, trial.d_high_load, *d_pparams);

         if (accept_break) {
            best_breakoff_box = corner_box;
            best_breakoff_size = corner_box_size;
            best_breakoff_load = corner_box_load;
         }

      } // while loop

      if (d_print_break_steps) {
         tbox::plog << "After corner #" << bn << ", best_breakoff_box is now box "
                    << bn << " " << best_breakoff_box
                    << ", best_breakoff_size = " << best_breakoff_size
                    << ", best_breakoff_load = " << best_breakoff_load
                    << std::endl;
      }

   } // bn loop

   if (!best_breakoff_box.empty()) {
      trial.computeBreakData(best_breakoff_box);
   }

#ifdef DEBUG_CHECK_ASSERTIONS
   for (hier::BoxContainer::iterator bi = trial.d_breakoff.begin();
        bi != trial.d_breakoff.end();
        ++bi) {
      const hier::Box& b = *bi;
      const hier::IntVector s = b.numberCells();
      for (int d = 0; d < dim.getValue(); ++d) {
         if (((s(d) < d_pparams->getMinBoxSize() (d)) && (s(d) != box_dims(d))) ||
             (s(d) > box_dims(d))) {
            TBOX_ERROR("BalanceBoxBreaker library error:\n"
               << "breakoff box " << b << ", with size " << s
               << "\nis not between the min size " << d_pparams->getMinBoxSize()
               << "\nand the original box size " << box_dims << "\n"
               << "orig box " << trial.d_whole_box << "\n"
               << "break box " << b << "\n"
               << "break box size " << b.size() << "\n"
               << "ideal brk load " << trial.d_ideal_load);
         }
      }
   }
   for (hier::BoxContainer::iterator bi = trial.d_leftover.begin();
        bi != trial.d_leftover.end();
        ++bi) {
      const hier::Box& b = *bi;
      const hier::IntVector s = b.numberCells();
      for (int d = 0; d < dim.getValue(); ++d) {
         if (((s(d) < d_pparams->getMinBoxSize() (d)) && (s(d) != box_dims(d))) ||
             (s(d) > box_dims(d))) {
            TBOX_ERROR("BalanceBoxBreaker library error:\n"
               << "leftover box " << b << ", with size " << s
               << "\nis not between the min size " << d_pparams->getMinBoxSize()
               << "\nand the original box size " << box_dims << "\n"
               << "orig box " << trial.d_whole_box << "\n"
               << "break box " << b << "\n"
               << "break box size " << b.size() << "\n"
               << "ideal brk load " << trial.d_ideal_load);
         }
      }
   }
#endif

   return !trial.d_breakoff.empty();
}

/*
 *************************************************************************
 *************************************************************************
 */

void
BalanceBoxBreaker::burstBox(
   hier::BoxContainer& boxes,
   const hier::Box& bursty,
   const hier::Box& solid)
{
   /*
    * This method lacks logic to handle the case of solid not being
    * completely inside bursty.  That feature is not currently needed.
    */
   TBOX_ASSERT(bursty.contains(solid));

   const hier::IntVector solid_size = solid.numberCells();

   boxes.clear();
   hier::Box cutme = bursty;
   while (!cutme.isSpatiallyEqual(solid)) {

      int cut_dir = 999999;
      bool cut_above_solid = false; // Whether to slice off the piece above solid (vs below).
      /*
       * Find direction and place to cut.  To minimize slivers, cut
       * from cutme the thickest slab (in direction normal to cut)
       * possible.
       */
      int slab_thickness = 0;
      for (tbox::Dimension::dir_t d = 0; d < solid_size.getDim().getValue(); ++d) {
         if (cutme.numberCells(d) > solid_size(d)) {
            const int thickness_from_upper_cut = cutme.upper() (d)
               - solid.upper() (d);
            if (thickness_from_upper_cut > slab_thickness) {
               slab_thickness = thickness_from_upper_cut;
               cut_dir = d;
               cut_above_solid = true;
            }
            const int thickness_from_lower_cut = solid.lower() (d)
               - cutme.lower() (d);
            if (thickness_from_lower_cut > slab_thickness) {
               slab_thickness = thickness_from_lower_cut;
               cut_dir = d;
               cut_above_solid = false;
            }
         }
      }
      TBOX_ASSERT(cut_dir >= 0 && cut_dir < solid_size.getDim().getValue());

      hier::Box removeme = cutme;
      if (cut_above_solid) {
         cutme.setUpper(static_cast<hier::Box::dir_t>(cut_dir),
            solid.upper(static_cast<hier::Box::dir_t>(cut_dir)));
         removeme.setLower(static_cast<hier::Box::dir_t>(cut_dir),
            solid.upper(static_cast<hier::Box::dir_t>(cut_dir)) + 1);
      } else {
         cutme.setLower(static_cast<hier::Box::dir_t>(cut_dir),
            solid.lower(static_cast<hier::Box::dir_t>(cut_dir)));
         removeme.setUpper(static_cast<hier::Box::dir_t>(cut_dir),
            solid.lower(static_cast<hier::Box::dir_t>(cut_dir)) - 1);
      }

      boxes.push_back(removeme);

   }
}

/*
 *************************************************************************
 * Compute a box_size score that is low for box widths smaller than
 * some threshold_width.  Boxes longer than this are okay assuming
 * that they can be broken down further.  Boxes shorter than this are
 * undesirable because they lead to excessive box count and/or aspect
 * ratios.
 *************************************************************************
 */
double BalanceBoxBreaker::computeWidthScore(
   const hier::IntVector& box_size,
   double threshold_width)
{
   return tbox::MathUtilities<double>::Min(box_size.min() / threshold_width, 1);
}

/*
 *************************************************************************
 * Compute a size score for multiple boxes.  The combined score is the
 * product of individual scores.
 *************************************************************************
 */
double BalanceBoxBreaker::computeWidthScore(
   const hier::BoxContainer& boxes,
   double threshold_width)
{
   double width_score = 1.0;
   for (hier::BoxContainer::const_iterator bi = boxes.begin(); bi != boxes.end(); ++bi) {
      width_score *= computeWidthScore(bi->numberCells(), threshold_width);
   }
   return width_score;
}

/*
 *************************************************************************
 *************************************************************************
 */
void BalanceBoxBreaker::setTimers()
{
   t_break_off_load = tbox::TimerManager::getManager()->
      getTimer("mesh::BalanceBoxBreaker::breakOffLoad()");
   t_find_bad_cuts = tbox::TimerManager::getManager()->
      getTimer("mesh::BalanceBoxBreaker::find_bad_cuts");
}

/*
 *************************************************************************
 * Set all members invariant with each trial break.
 *************************************************************************
 */
BalanceBoxBreaker::TrialBreak::TrialBreak(
   const PartitioningParams& pparams,
   double threshold_width,
   const hier::Box& whole_box,
   double whole_box_load,
   const std::vector<double>& corner_weights,
   const std::vector<std::vector<bool> >& bad_cuts,
   double ideal_load,
   double low_load,
   double high_load):
   d_breakoff_load(0.0),
   d_breakoff(),
   d_leftover(),
   d_ideal_load(ideal_load),
   d_low_load(low_load),
   d_high_load(high_load),
   d_width_score(1.0),
   d_balance_penalty(0.0),
   d_pparams(&pparams),
   d_threshold_width(threshold_width),
   d_whole_box(whole_box),
   d_whole_box_load(whole_box_load),
   d_corner_weights(corner_weights),
   d_bad_cuts(bad_cuts)
{
   computeMerits(); // This is the merits of doing nothing.
}

/*
 *************************************************************************
 *************************************************************************
 */
BalanceBoxBreaker::TrialBreak::TrialBreak(
   const TrialBreak& orig,
   bool make_reverse):
   d_breakoff_load(0.0),
   d_breakoff(),
   d_leftover(),
   d_ideal_load(orig.d_whole_box_load - orig.d_ideal_load),
   d_low_load(orig.d_whole_box_load - orig.d_high_load),
   d_high_load(orig.d_whole_box_load - orig.d_low_load),
   d_width_score(orig.d_width_score),
   d_balance_penalty(orig.d_balance_penalty),
   d_pparams(orig.d_pparams),
   d_threshold_width(orig.d_threshold_width),
   d_whole_box(orig.d_whole_box),
   d_whole_box_load(orig.d_whole_box_load),
   d_corner_weights(orig.d_corner_weights),
   d_bad_cuts(orig.d_bad_cuts)
{
   // make_reverse only prevents this from being mistaken for a copy constructor.
   NULL_USE(make_reverse);
   computeMerits(); // This is the merits of doing nothing.
}

/*
 *************************************************************************
 * Break box out of d_whole_box and store results.
 *************************************************************************
 */
void BalanceBoxBreaker::TrialBreak::computeBreakData(
   const hier::Box& box)
{
   d_breakoff.clear();
   d_leftover.clear();
   burstBox(d_leftover, d_whole_box, box);
   bool do_cut = box.size() >= d_pparams->getMinimumCellRequest();
   if (do_cut) {
      for (auto itr = d_leftover.begin(); itr != d_leftover.end(); ++itr) {
         if ((*itr).size() < d_pparams->getMinimumCellRequest()) {
            do_cut = false;
         }
      }
   }
   if (do_cut) {
      d_breakoff.push_back(box);
      d_breakoff_load = computeBreakOffLoad(box);
   } else {
      d_leftover.clear(); 
      d_leftover.push_back(d_whole_box);
      d_breakoff_load = 0;
   }
   computeMerits();
}

/*
 *************************************************************************
 * Compute the load total represent by the Box if it is broken off of
 * d_whole_box.
 *************************************************************************
 */
double BalanceBoxBreaker::TrialBreak::computeBreakOffLoad(
   const hier::Box& box)
{
   double breakoff_load = 0.0;

   if (d_corner_weights.empty()) {
      /*
       * If there is no corner weight information, the breakoff load is
       * box.size() multiplied by the ratio of the whole load to the whole
       * box size.
       */
      breakoff_load = (d_whole_box_load /
                       static_cast<double>(d_whole_box.size())) *
                      static_cast<double>(box.size());
   } else {
      /*
       * If their are corner weights, use them in the load computation.
       */
      const tbox::Dimension& dim = box.getDim();
      hier::IntVector mid_point(d_whole_box.numberCells()/2);
      hier::IntVector corner_id(hier::IntVector::getZero(dim));
      unsigned int i = 0;
      /*
       * loop over all the corners.
       */
      do {
         /*
          * d_corner_box stores the corner boxes for repeated calls of this
          * method.
          */
         if (d_corner_box.size() == i) {
            d_corner_box.push_back(d_whole_box);
            hier::Box& current_box = d_corner_box[i];
            for (unsigned int d = 0; d < dim.getValue(); ++d) {
               if (corner_id[d] == 0) {
                  current_box.setLower(d, d_whole_box.lower()[d]);
                  current_box.setUpper(d, d_whole_box.lower()[d] + mid_point[d] - 1);
               } else {
                  current_box.setLower(d, d_whole_box.lower()[d] + mid_point[d]);
                  current_box.setUpper(d, d_whole_box.upper()[d]);
               }
            }
         }

         /*
          * Increment n-dimensional corner id
          */
         for (unsigned int d = 0; d < dim.getValue(); ++d) {
            if (corner_id[d] == 0) {
               corner_id[d] = 1;
               break;
            } else {
               corner_id[d] = 0;
            }
         }

         /*
          * Compute the load for the intersection of the input box and
          * the current corner box.  If the intersection is equal to the
          * input box, we don't need to look at any more corner boxes and
          * can break out of the do-while loop.
          */
         hier::Box intersection(box * d_corner_box[i]);
         if (!intersection.empty()) {
            double corner_load = d_corner_weights[i] * d_whole_box_load;
            breakoff_load += (corner_load /
                              static_cast<double>(d_corner_box[i].size())) *
                             static_cast<double>(intersection.size());
            if (intersection.isSpatiallyEqual(box)) {
               break;
            }
         }
         ++i;
      } while (corner_id != hier::IntVector::getZero(dim));
   }
   return breakoff_load;
}


/*
 *************************************************************************
 *************************************************************************
 */
void BalanceBoxBreaker::TrialBreak::swapWithReversedTrial(
   TrialBreak& reversed)
{
#ifdef DEBUG_CHECK_ASSERTIONS
   if (&d_whole_box != &reversed.d_whole_box ||
       &d_bad_cuts != &reversed.d_bad_cuts ||
       d_pparams != reversed.d_pparams ||
       d_threshold_width != reversed.d_threshold_width) {
      TBOX_ERROR("BalanceBoxBreaker::TrialBreak::swapWithReversedTrial:"
         << "\nIncompatible TrialBreaks:"
         << "\nthis:\n" << *this
         << "\nreversed:\n" << reversed);
   }
#endif

   d_breakoff.swap(reversed.d_leftover);
   d_leftover.swap(reversed.d_breakoff);
   d_corner_box.swap(reversed.d_corner_box);
   d_breakoff_load = d_whole_box_load - reversed.d_breakoff_load;
   computeMerits();
}

/*
 *************************************************************************
 *************************************************************************
 */
bool BalanceBoxBreaker::TrialBreak::computeMerits()
{
   d_width_score =
      (d_breakoff.empty() && d_leftover.empty()) ?
      computeWidthScore(d_whole_box.numberCells(), d_threshold_width) :
      computeWidthScore(d_breakoff, d_threshold_width)
      * computeWidthScore(d_leftover, d_threshold_width);

   d_balance_penalty = BalanceBoxBreaker::computeBalancePenalty(
         d_breakoff_load - d_ideal_load);

   return BalanceUtilities::compareLoads(
      d_flags, 0, d_breakoff_load,
      d_ideal_load, d_low_load, d_high_load, *d_pparams);
}

/*
 *************************************************************************
 * Determine if this TrialBreak is an improvement over another.  If
 * both are in-range, pick the one with the best width score.  If only
 * one is in-range, pick that one.  If none are in range but both are
 * better than not breaking, chose the one with better load.  If only
 * one is better than not breaking, choose that one.  Else, choose
 * none.
 *************************************************************************
 */
int BalanceBoxBreaker::TrialBreak::improvesOver(
   const TrialBreak& other) const
{
   int improves = 0;

   if (this->d_flags[3] && other.d_flags[3]) {
      improves =
         (this->d_width_score > other.d_width_score) ? 1 :
         (this->d_width_score < other.d_width_score) ? -1 : 0;
      if (improves == 0) {
         improves =
            (this->d_balance_penalty < other.d_balance_penalty) ? 1 :
            (this->d_balance_penalty > other.d_balance_penalty) ? -1 : 0;
         // How about balance_score = -balance_penalty so everything is a "score."
      }
   } else if (this->d_flags[3]) {
      improves = 1;
   } else if (other.d_flags[3]) {
      improves = -1;
   } else if (this->d_flags[2] == 1 && other.d_flags[2] == 1) {
      int flags[4] = { 0, 0, 0, 0 };
      improves = BalanceUtilities::compareLoads(
            flags, other.d_breakoff_load, this->d_breakoff_load,
            d_ideal_load, d_low_load, d_high_load, *d_pparams);
   } else if (this->d_flags[2] == 1) {
      improves = 1;
   } else if (other.d_flags[2] == 1) {
      improves = -1;
   }

   return improves == 1;
}

/*
 *************************************************************************
 *************************************************************************
 */
void BalanceBoxBreaker::TrialBreak::swap(TrialBreak& other)
{
   TBOX_ASSERT(&d_whole_box == &other.d_whole_box);
   TBOX_ASSERT(&d_bad_cuts == &other.d_bad_cuts);
   TBOX_ASSERT(d_pparams == other.d_pparams);
   TBOX_ASSERT(d_ideal_load == other.d_ideal_load);
   TBOX_ASSERT(d_low_load == other.d_low_load);
   TBOX_ASSERT(d_high_load == other.d_high_load);
   TBOX_ASSERT(d_threshold_width == other.d_threshold_width);

   d_breakoff.swap(other.d_breakoff);
   d_leftover.swap(other.d_leftover);
   d_corner_box.swap(other.d_corner_box);

   double tmpd;

   tmpd = d_breakoff_load;
   d_breakoff_load = other.d_breakoff_load;
   other.d_breakoff_load = tmpd;

   tmpd = d_width_score;
   d_width_score = other.d_width_score;
   other.d_width_score = tmpd;

   tmpd = d_balance_penalty;
   d_balance_penalty = other.d_balance_penalty;
   other.d_balance_penalty = tmpd;

   int tmpi;
   for (int i = 0; i < 4; ++i) {
      tmpi = d_flags[i];
      d_flags[i] = other.d_flags[i];
      other.d_flags[i] = tmpi;
   }
}

/*
 *************************************************************************
 *************************************************************************
 */
std::ostream& operator << (
   std::ostream& co,
   const BalanceBoxBreaker::TrialBreak& tb)
{
   co.unsetf(std::ios::fixed | std::ios::scientific);
   co.precision(6);
   co << "\n      TrialBreak state: broke off "
   << tb.d_breakoff_load << " / " << tb.d_ideal_load
   << " [" << tb.d_low_load << ", " << tb.d_high_load
   << "] from " << tb.d_whole_box << '|'
   << tb.d_whole_box.numberCells() << '|'
   << tb.d_whole_box.size() << " into:\n"
   << "      " << " breakoff boxes: " << tb.d_breakoff.format("           ")
   << "      " << " leftover boxes: " << tb.d_leftover.format("           ")
   << "      imbalance: "
   << (tb.d_breakoff_load - tb.d_ideal_load)
   << " d_balance_penalty: " << tb.d_balance_penalty
   << " d_width_score: " << tb.d_width_score
   << "\n      d_flags:"
   << "  " << tb.d_flags[0]
   << "  " << tb.d_flags[1]
   << "  " << tb.d_flags[2]
   << "  " << tb.d_flags[3]
   << std::endl;
   return co;
}

}
}

#if !defined(__BGL_FAMILY__) && defined(__xlC__)
/*
 * Suppress XLC warnings
 */
#pragma report(enable, CPPC5334)
#pragma report(enable, CPPC5328)
#endif

#endif
